A Levinson-like algorithm for symmetric strongly nonsingular higher order semiseparable plus band matrices

In this paper we will derive a solver for a symmetric strongly nonsingular higher order generator representable semiseparable plus band matrix. The solver we will derive is based on the Levinson algorithm, which is used for solving strongly nonsingular Toeplitz systems. In a first part an O(p 2 n) solver for a semiseparable matrix of semiseparability rank p is derived, and in a second part we derive an O(l 2 n) solver for a band matrix with bandwidth 2l + 1. Both solvers are constructed in a similar way: firstly a Yule-Walker-like equation needs to be solved, and secondly this solution is used for solving a linear equation with an arbitrary right-hand side. Finally a combination of the above methods is presented to solve linear systems with semiseparable plus band coefficient matrices. The overall complexity of this solver is 6(l + p) 2 n plus lower order terms. In a final section numerical experiments are performed. Attention is paid to the timing and the accuracy of the described methods. A Levinson-like algorithm for symmetric strongly nonsingular higher order semiseparable plus band matrices * Abstract In this paper we will derive a solver for a symmetric strongly nonsingular higher order generator rep-resentable semiseparable plus band matrix. The solver we will derive is based on the Levinson algorithm, which is used for solving strongly nonsingular Toeplitz systems. In a first part an O(p 2 n) solver for a semiseparable matrix of semiseparability rank p is derived, and in a second part we derive an O(l 2 n) solver for a band matrix with bandwidth 2l + 1. Both solvers are constructed in a similar way: firstly a Yule-Walker-like equation needs to be solved, and secondly this solution is used for solving a linear equation with an arbitrary right-hand side. Finally a combination of the above methods is presented to solve linear systems with semiseparable plus band coefficient matrices. The overall complexity of this solver is 6(l + p) 2 n plus lower order terms. In a final section numerical experiments are performed. Attention is paid to the timing and the accuracy of the described methods.